/*
 * nqueens.c
 *
 * Author: Michael Stegeman
 * Date: 4/17/08
 *
 * Description: Program to solve N-queens problem using simulated annealing
 *
 * Usage: Change the #define N to whatever number you want to solve, make, and run
 *
 *
 *
 * The MIT License
 *
 * Copyright (c) <2008> <Michael Stegeman>
 *
 * Permission is hereby granted, free of charge, to any person obtaining a copy
 * of this software and associated documentation files (the "Software"), to deal
 * in the Software without restriction, including without limitation the rights
 * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
 * copies of the Software, and to permit persons to whom the Software is
 * furnished to do so, subject to the following conditions:
 *
 * The above copyright notice and this permission notice shall be included in
 * all copies or substantial portions of the Software.
 *
 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
 * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
 * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
 * THE SOFTWARE.
 *
 */

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <math.h>
#include <string.h>
#define N 20 /* must be >=4 */

/* Searches array s for value row, return 1 if found, else return 0 */
int find(int *s, int cur, int row)
{
    int i;

    for(i = 0; i < cur; ++i)
    {
        if(s[i] == row)
        {
            return 1;
        }
    }

    return 0;
}

/* Creates random initial state */
void initial_state(int *s)
{
    int i, row;
    
    for(i = 0; i < N; ++i)
    {
        do /* find a random row for each queen */
        {
            row = rand() % N;
        } while(find(s, i, row) == 1);

        s[i] = row;
    }
}

/* Makes a small, random change to current state */
int * perturb(int *s)
{
    int *s_nearby;
    int i, j, temp;
    int iterations = 0;

    while(1)
    {
        /* pick 2 random queens to switch */
        i = rand() % N;
        j = rand() % N;

        /* Check if new state is immediately invalid, i.e. an immediate diagonal 
         * puts the queen in danger, makes for better switches, or that random ints are equal.
         * Also, limit the number of iterations, to avoid getting stuck in a random number loop */
        if((iterations != 10) && ((i == j) || \
            ((i == 0) && ((s[i+1] == (s[j] + 1)) || (s[i+1] == (s[j] - 1)))) || \
            ((i == (N - 1)) && ((s[i-1] == (s[j] + 1)) || (s[i-1] == (s[j] - 1)))) || \
            ((j == 0) && ((s[j+1] == (s[i] + 1)) || (s[j+1] == (s[i] - 1)))) || \
            ((j == (N - 1)) && ((s[j-1] == (s[i] + 1)) || (s[j-1] == (s[i] - 1)))) || \
            (((i > 0) && (i < (N -1))) && ((s[i+1] == (s[j] + 1)) || (s[i+1] == (s[j] - 1)) || (s[i-1] == (s[j] + 1)) || (s[i-1] == (s[j] - 1)))) || \
            (((j > 0) && (j < (N -1))) && ((s[j+1] == (s[i] + 1)) || (s[j+1] == (s[i] - 1)) || (s[j-1] == (s[i] + 1)) || (s[j-1] == (s[i] - 1))))))
        {
            ++iterations;
            continue;
        }
        else
        {
            break;
        }
    }

    /* switch queens */
    s_nearby = (int *)malloc(N * sizeof(int));
    memcpy(s_nearby, s, N * sizeof(int)); /* copy original array into new array */
    temp = s_nearby[i];
    s_nearby[i] = s_nearby[j];
    s_nearby[j] = temp;

    return s_nearby;
}

/* Calculates "energy" of current state, aka number of improperly placed queens */
int E(int *s)
{
    int bad = 0;
    int i, j, k;

    for(i = 0; i < (N - 1); ++i) /* loop through columns of chess board */
    {
        j = i + 1;
        for(k = s[i] - 1; k > -1; --k) /* loop through diagonals forward and upward */
        {
            if(j >= N)
            {
                break;
            }
            else if(s[j] == k) /* check if a queen occupies current space */
            {
                ++bad;
            }

            ++j;
        }

        j = i + 1;
        for(k = s[i] + 1; k < N; ++k) /* loop through diagonals forward and downward */
        {
            if(j >= N)
            {
                break;
            }
            else if(s[j] == k) /* check if a queen occupies current space */
            {
                ++bad;
            }

            ++j;
        }
    }

    return bad; /* return number of bad queens */
}

/* Calculates the probability that the current state should change
 * based on current energy, energy of nearby state, and temperature */
float P(int e, int e_nearby, float T)
{
    int delta_e = abs(e_nearby - e); /* find change in energy */
    float p = expf(-(delta_e / T)); /* calculate probability that we should change */

    /* if e_nearby is lower than the current energy, slightly raise the probability */
    if(e_nearby < e)
    {
        p += 0.1;
    }
    
    return p;
}

/* Prints N-queens solution */
void print_board(int *s)
{
    int i, j;

    fprintf(stdout, "Solution for %dx%d Puzzle:\n", N, N);
    for(i = 0; i < N; ++i) /* loop through row by row */
    {
        for(j = 0; j < N; ++j) /* loop through column by column */
        {
            if(s[j] == i) /* print a Q if there is a queen in the space */
                fprintf(stdout, "Q ");
            else
                fprintf(stdout, "~ ");
        }

        printf("\n");
    }
}

int main()
{
    float T = 5 * N; /* Starting temperature */
    float Tmax = T;
    int *s = (int *)malloc(sizeof(int) * N);
    int *s_nearby;
    int e, e_nearby;

    srand((unsigned)time(NULL)); /* seed rand() function */
    initial_state(s); /* start with some random initial state */
    e = E(s); /* e is energy of current state */

    /* Loop until temperature is near freezing */
    while (T > 0)
    {
        /* examine a nearby state */
        s_nearby = perturb(s);
        e_nearby = E(s_nearby);

        /* check for 0 energy (found a solution) */
        if(e_nearby == 0)
        {
            free(s);
            s = s_nearby;
            s_nearby = NULL;
            break;
        }
        
        // decide whether to move to new state
        if (P(e, e_nearby, T) > ((float)rand() / (float)RAND_MAX))
        {
            free(s);
            s = s_nearby;
            s_nearby = NULL;
            e = e_nearby;
        }
        else
        {
            free(s_nearby);
            s_nearby = NULL;
        }

        T = T * 0.999999; /* temperature drops a bit */

        /* raise temperature to avoid sitting in a local minimum energy state */
        if(T < 0.1)
        {
            T = Tmax * 0.6667;
            Tmax = T;
        }
    }

    /* solution found, print it */
    print_board(s);

    /* free memory */
    free(s);
    s = NULL;
}

